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NUMERICAL  SIMULATION  OF  TWO-DIMENSIONAL 
SPATIALLY-DEVELOPING  MIXING  LAYERS 

R.  V.  Wilson^  and  A.  O.  Demuren 
Department  of  Mechanical  Engineering 
Old  Dominion  University 
Norfolk.  VA  23529 


ABSTRACT 


Two-dimensional,  incompressible,  spatially  developing  mixing  layer  simulations  are  per¬ 
formed  at  Re  =  10^  and  1(1*  with  two  classes  of  perturbations  applied  at  the  inlet  boundary;  (0 
combinations  of  discrete  modes  from  linear  stability  theory,  and  (I’O  a  broad  spectrum  of  modes 
derived  from  experimentally  measured  velocity  spectra.  The  effect  of  the  type  and  strength  of 
inlet  perturbations  on  vortex  dynamics  and  time-averaged  properties  are  explored.  Two-point  spa¬ 
tial  velocity  and  autocorrelations  are  used  to  estimate  the  size  and  lifetime  of  the  resulting  coher¬ 
ent  structures  and  to  explore  possible  feedback  effects.  The  computed  time-averaged  properties 
such  as  mean  velocity  probles,  turbulent  statistics,  and  spread  rates  show  good  agreement  with 
experimentally  measured  values.  It  is  shown  that  by  forcing  with  a  broad  spectrum  of  modes 
derived  from  an  experimental  energy  spectrum  many  experimentally  observed  phenomena  can  be 
reproduced  by  a  2-D  simulation.  The  strength  of  the  forcing  merely  affected  the  length  required 


for  the  dominant  coherent  structures  to  become  fully-developed.  Thus  intensities  comparable  to 


those  of  the  background  turbulence  in  many  wind 
given  sufficient  simulation  length. 


tunnel  experiments  produced  the  same  results,  3^ 
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L  INTRODUCTION 


Mixing  layers  are  studied  extensively  because  of  their  wide  practical  applications  as  well  as 
for  scientific  interest.  Apart  from  the  obvious  cases  of  two  stream  mixing  layers  and  the  near  field 
of  jets  they  can  be  found  in  most  shear  flows  away  from  walls.  Many  complicated  flow  situations 
can  be  shown  to  be  a  synthesis  of  wall  boundary  layers  and  mixing  layers,  amongst  others.  Thus, 
even  in  recirculating  flows,  mixing  layers  play  an  important  role  and  their  understanding  is  neces¬ 
sary  for  the  prediction  of  many  practical  flow  situations.  Of  the  two  major  flow  types:  the  wall 
boundary  layer  and  the  free  mixing  layer,  the  latter  is  the  less  understood,  more  controversial  and 
less  predictable. 

A  few  of  the  major  controversies  are  the  importance  of  coherent  structures  in  free  shear  flows, 
their  dimensionality  and  their  relationship  to  real  turbulence.  Following  Brown  and  Roshko^’s 
landmark  study  on  the  dominant  role  of  coherent  structures  in  the  evolution  of  mixing  layers,  and 
their  persistence  even  at  high  Reynolds  numbra^s,  several  experimental  investigations^"*  have 
been  performed  to  explore  the  role  of  the  coherent  structures,  the  mechanisms  for  their  formation 
and  their  eventual  breakdown  into  random  turbulent  motion.  Winant  and  Browand^  observed  that 
the  process  of  vortex  pairing  was  responsible  for  shear  layer  growth.  This  was  essentially  a  two- 
dimensional  process.  Dimotakis  and  Brown^  confirmed  the  persistence  of  the  coherent  structures 
at  high  Reynolds  numbers  and  that  large-scale  entrainment  and  small-scale  mixing  occurred 
through  distinct  processes  at  different  stages  of  evolution  of  the  coherent  structures.  Ho  and 
Huang^  studied  mixing  layers  perturbed  with  disturbances  of  different  frequencies.  They  found 
that  vortex  pairing  was  controlled  by  the  presence  of  subharmonic  modes  of  the  fundamental  fre¬ 
quency  of  the  mixing  layer.  Thus  vortex  mergings  at  specific  locations  involving  a  desired  number 
of  vortices  could  be  initiated  simply  by  forcing  the  layer  with  disturbances  which  contain  corre- 
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spending  subharmonic  modes.  Later  results  by  Ho  et  al.^  have  downplayed  the  ability  to  control 
this  process  because  phase  jitter  could  lead  to  decorrelation  between  the  fundamental  and  its  sub¬ 
harmonic  modes  leading  to  more  random  merging  locations.  Perhaps  the  major  challenge  to  the 
Brown  and  Roshko  picture  of  the  role  of  coherent  struemres  was  presented  by  Chandrsuda  et  al.^. 
Based  on  the  analysis  of  their  experimental  data,  they  postulated  that  this  must  be  a  relic  of  the 
transition  process  and  that  the  spanwise  vortices  would  only  be  dominant  in  the  absence  of  ini¬ 
tially  turbulent  boundary  layers  at  the  splitter  plate,  high  free  stream  turbulence,  or  finite  three- 
dimensional  disturbances.  This  assertion  was  partly  confirmed  by  Browand  and  Latigo^  in  their 
study  of  a  plane  mixing  layer,  but  other  studies  such  as  that  by  Hussain  and  Zaman*  contradict 
this  finding. 

It  is  now  believed^  that  a  feedback  mechanism  could  produce  a  disturbance,  with  a  frequency 
related  to  the  size  of  the  domain,  which  would  interact  with  the  mixing  layer  to  influence  the 
development  all  the  way  to  the  inlet.  Thus,  the  mixing  layer  could  never  forget  its  origin  and  self¬ 
similarity  conditions  would  not  exist.  In  a  comprehensive  review  of  experimental  data  of  single 
stream  mixing  layers  developing  from  different  initial  conditions,  Birch^  found  that  for  turbulent 
initial  conditions  the  spread  rate  in  the  developing  region  could  overshoot  the  asymptotic  value  in 
the  fully-developed  region  by  over  20%,  whereas  with  laminar  initial  conditions  the  spread  rate 
monotonically  increased  to  its  asymptotic  value  with  no  associated  overshoot.  The  author  also 
found  that  for  turbulent  initial  conditions,  the  turbulent  shearing  stress  gradually  approached  the 
far-field,  asymptotic  value,  while  with  laminar  initial  conditions  an  overshoot  twice  the  asymp¬ 
totic  value  of  the  turbulent  shearing  stress  occurred.  Bradshaw^®  reported  the  persistence  of  initial 
conditions  up  to  1000  initial  momentum  thicknesses  in  the  mixing  region  of  a  round  jet 
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Linearized  stability  theory  has  been  successful  in  predicting  the  initial  development  of  the 
shear  layer  before  the  onset  of  vortex  interactions  which  are  in  general  strongly  nonlinear. 
Michalke^^  performed  a  temporal  stability  analysis  of  a  single  plane  mixing  layer,  which  showed 
a  finite  band  of  wavenumbers  experiencing  exponential  growth  rates.  All  modes  experienced  con¬ 
stant  phase  speed  independent  of  wavenumber.  A  stability  analysis  of  the  spatially-developing 
plane  mixing  layer  was  performed  by  Michalke^^,  and  Monkewitz  and  Huerre*^,  with  the  latter 
considering  both  the  Blasius  and  hyperbolic  tangent  velocity  profiles  of  arbitrary  velocity  ratio. 
For  both  velocity  profiles  it  was  found  that  the  maximum  spatial  growth  rate  was  approximately 
proportional  to  the  velocity  ratio,  with  a  finite  range  of  frequencies  that  experienced  exponential 
growth  rates.  The  phase  velocity  was  found  to  be  a  function  of  frequency,  with  frequencies  less 
than  the  most  amplified  mode  experiencing  a  larger  phase  speed.  The  phase  speed  was  approxi¬ 
mately  constant  and  equal  to  that  predicted  by  temporal  analysis  for  modes  greater  than  the  most 
amplified  mode.  It  was  found  that  linearized  stability  theory  applied  to  a  spatially  developing 
mixing  layer  yielded  growth  rates  which  showed  better  agreement  with  experimental  values. 

Linear  theory  is  not  capable  of  predicting  strongly  nonlinear  behavior  such  as  saturation  of  the 
growth  rates,  suppression  of  one  mode  by  another,  and  interaction  of  two  or  more  vortices.  The 
above  mentioned  phenomena  must  be  simulated  numerically.  Riley  and  Metcalfe*"^  performed  2- 
D  and  3-D  temporally-developing  mixing  layer  simulations  and  found  that  when  forced  with  a 
broad  range  of  modes  plus  the  fundamental  mode  at  strong  amplitudes,  the  shear  layer  exhibited 
larger  initial  growth  rates  followed  by  suppressed  growth  rates  compared  to  simulations  forced 
with  only  a  broad  range  of  modes.  Cain  et  al.^^  studied  the  temporally-developing  mixing  layer 
and  noted  that  by  forcing  the  initial  conditions  with  fluctuations  at  relatively  high  intensity,  the 
initial  growth  rate  was  two  to  three  times  the  asymptotic  value.  It  was  also  found  that  the  growth 
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of  the  momentum  thickness  was  accelerated  by  concentrating  the  energy  of  the  initial  conditions 
in  the  most  amplified  band  of  wavenumbers  predicted  by  linear  theory.  Corcos  and  Sherman*^ 
performed  two-dimensional  (2-D),  time-developing  mixing  layer  simulations  and  reported  that 
the  resulting  vortex  dynamics  were  largely  determined  by  the  subharmonic  frequency  content  in 
the  initial  conditions.  With  only  the  fundamental  mode  present  in  the  initial  conditions,  the  ini¬ 
tially  uniform  vortex  sheet  of  the  shear  layer  rolled  up  into  an  array  of  spanwise  vorticity  “lumps” 
separated  by  a  distance  of  one  fundamental  wavelength.  The  addition  of  a  subharmonic  mode  in 
phase  with  that  of  the  fundamental  mode  resulted  in  two  spanwise  vorticity  lumps  merging  to 
form  a  single  larger  structure  which  caused  a  sudden  growth  of  the  shear  layer  width.  Corcos  and 
Lm  performed  linearized  3-D  simulations  and  reported  a  hierarchy  of  motions  in  which  the  2-D 
roll-up  and  pairing  of  spanwise  vortices  drove  the  secondary  instabilities  to  produce  streamwise 
vortices  in  the  braid  region  of  the  2-D  structures.  Three-dimensional  effects  are  further  discussed 
by  Metcalfe  et  al.'*,  Lesieur  et  al.^^,  and  Rogers  and  Moser^®. 

Davis  and  Moore^*  studied  the  manipulation  of  vortex  pairing  through  the  forcing  of  2-D  axi- 
symmetric  and  plane  spatially-developing  shear  layers.  The  authors  showed  that  the  vortex  pair¬ 
ing  events  were  determined  by  inlet  boundary  conditions  with  the  number  of  pairings  equal  to  the 
number  of  subharmonic  modes  present  in  the  inlet  boundary  conditions.  The  vortex  dynamics 
were  found  to  be  relatively  independent  of  Reynolds  number.  Lowery  and  Reynolds^^  examined 
the  effect  of  vortex  pairing  on  time-averaged  properties  by  performing  highly  resolved  2-D  and  3- 
D  spatial  DNS  of  the  forced  mixing  layer.  The  forcing  function  for  the  2-D  simulations  consisted 
of  the  fundamental,  and  first  and  second  subharmonic  modes  with  no  explicit  phase  shift.  Around 
each  region  of  vortex  pairing,  the  layer  experienced  rapid  growth,  and  positive  turbulence  kinetic 
energy  production,  while  between  pairing  events,  inhibited  layer  growth  and  negative  tiu'bulent 
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kinetic  energy  production  occiured.  The  2-D  simulations  produced  spread  rates  consistent  with 
those  in  the  3-D  simulations  and  experiment. 

The  literature  shows  that  temporal  simulations  are  capable  of  predicting  experimentally 
observed  phenomena  such  as  shear  layer  rollup,  vortex  pairing,  as  well  as  experimentally  mea¬ 
sured  spread  rates.  However,  the  spatial  reference  frame  must  be  adopted  to  realize  asymmetric 
entrainment  in  which  a  greater  proportion  of  high  speed  fluid  is  entrained  into  the  layer,  disper¬ 
sion  of  instability  waves  due  to  mode  dependent  phase  speed  as  predicted  by  spatial  stability  the¬ 
ory,  and  feedback  effects  whereby  downstream  events  influence  the  development  of  upstream 
events. 

In  the  present  study,  effects  of  Reynolds  number  and  inlet  forcing  functions  are  investigated 
with  a  2-D  simulation  of  the  spatially  evolving  mixing  layer.  Several  inlet  conditions  are  consid¬ 
ered,  based  on  discrete  modes  from  linear  stability  theory  and  based  on  energy  spectra  and  inten¬ 
sity  from  experimental  observations.  Combinations  of  both  discrete  and  broad  band  mode  forcing 
are  also  utilized  to  investigate  the  complex  interactions  which  occur  in  experimental  attempts  to 
control  shear  layer  growth  through  large  amplitude  forcing  of  turbulent  jets.  A  unified  treatment  is 
presented  of  which  can  be  used  to  explain  many  phenomena  which  have  been  observed  in  several 
different  studies,  some  apparently  in  conflict.  The  need  for  some  forcing  of  the  inlet,  and  the  type 
required  in  order  to  reproduce  experimental  data  are  discussed.  In  addition,  two-point  spatial  and 
auto-correlations  are  used  to  determine  the  size  of  the  coherent  structures. 
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n.  GOVERNING  EQUATIONS 


The  equations  of  motion  governing  two-dimensional,  incompressible,  unsteady  fluid  flow  are: 


du- 

dXj 


,  a  ,  ,  idp  aw, 

+3— (U,«;)  =  --^+V3— 5— 

dt  aXj  1  pdx,  dxjdxj 


(1) 


(2) 


where  the  unknown  variables  are  w  =  u/,  v  =  «2>  and  p,  which  are  the  x-  and  y-component  veloc¬ 
ities,  and  pressure,  respectively.  No  spatial  filtering  of  the  governing  equations  is  performed 
because  the  flow  is  dominated  by  the  large  scales.  The  assumption  of  two-dimensionality  impedes 
the  formation  of  small  scales  via  vortex  stretching,  whereas  the  smallest  scales  which  may  be 
present  in  the  inlet  conditions  are  dissipated.  The  coordinate  system,  flow  parameters,  and  compu¬ 
tational  domain  used  in  the  simulations  are  shown  in  Fig.  1. 


Simulations  were  started  with  initial  conditions  consisting  of  the  streamwise  component  of 
velocity,  m,  set  to  the  same  reference  hyperbolic  tangent  profile  everywhere  in  the  computational 
domain,  and  v  =  p  =  0.  The  domain  of  interest  for  the  present  study  starts  at  some  distance  down¬ 
stream  of  the  splitter  plate  where  the  high  and  low  speed  streams  have  merged  to  form  a  mono¬ 
tonic  velocity  profile.  The  hyperbolic  tangent  profile  is  adopted  in  the  current  study,  which  defines 
the  mean  streamwise  velocity  component,  C//,  at  the  inlet  as: 


U^(y)  =  Uc{l  + 


(3) 


where  y*  =  y  -  H/2,  Uc  =  5  is  ih®  convective  velocity,  X  =  (!/«  - 1//,)  /  (!/«  +  Ui)  is  the 

velocity  ratio,  and  8,^^  =  (du/dy)  is  the  vorticity  thickness  at  the  inlet  cross  section 

(denoted  by  subscript  “o”),  with  (du/dy)  being  the  maximum  derivative  at  the  inlet  cross-sec- 


tion.  The  vorticity  thickness  is  a  measure  of  the  shear  layer  width  and  for  this  study,  =  0.1359 
is  used.  For  the  hyperbolic  tangent  profile,  the  initial  momentum  thickness  is  given  by, 
=  5h,/4  .  The  mean  flow  at  the  inlet  is  approximated  as  parallel,  resulting  in  a  zero  mean  trans¬ 
verse  velocity,  (y)  =0.  The  zero-gradient  boundary  condition  is  used  in  all  cases  for  the  speci¬ 
fication  of  the  instantaneous  pressure  at  the  inlet  boundary,  dp/dx,  =  0. 

Two  classes  of  perturbations  are  used  in  the  current  study:  (j)  perturbations  derived  from  lin¬ 
earized  inviscid  stability  analysis,  and  (ii)  perturbations  derived  from  experimentally  measured 
velocity  spectra  and  transverse  root  mean  square  (ms)  values.  The  inlet  boundary  condition  is 
then  a  superposition  of  fl\e  mean  and  perturbed  velocity  fields.  The  two  classes  of  perturbations 
are  defined  below. 

(i)  Perturbations  from  linearized  inviscid  stability  analysis 

The  first  class  of  perturbations  is  derived  from  the  solution  of  the  Rayleigh  equation  govern¬ 
ing  the  instability  of  the  reference  hyperbolic  tangent  profile  to  spatially  developing  disturbances. 
The  solution  yields  the  set  of  complex  eigenfunctions  corresponding  to  a  specified  frequency,  . 
In  this  study,  various  combinations  of  the  fundamental  (the  mode  which  invisi  id  theory  predicts 
the  most  amplified  growth  rate),  and  the  first  and  ^ond  subharmonic  modes  are  used  as  pertur¬ 
bations.  In  simulations  forced  with  more  than  one  discrete  mode,  no  explicit  phase  shift  was  intro¬ 
duced  between  the  modes.  A  complete  derivation  and  solution  of  the  Rayleigh  equation  is 
provided  by  Monkewitz  and  Huene^^,  with  the  general  form  of  the  perturbation  being; 

«V(y,0  =  Real{Kf(y)e‘V}  (4) 
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with  i  =  ,  t  being  the  time,  and  uf  (y)  being  the  complex  eigenfunction  for  the  mode  with 

angular  frequency,  The  eigenfunctions  for  the  fundamental,  first,  and  second  subharmonic 
modes,  were  obtained  from  Lowery  and  Reynolds^^. 

(ii)  Random  perturbations  from  experimental  data 

Perturbations  having  a  broad  spectrum  resembling  that  of  fully-developed,  mostly  random,  3- 
D  turbulence  were  generated.  The  perturbations  are  typical  of  those  found  in  the  experimental 
mixing  layer  originating  from  turbulent  boundary  layers  upstream  of  the  splitter  plate,  with  the 
obvious  difference  that  the  conditions  generated  here  are  2-D,  while  those  of  the  experiment  are 
3-D,  The  velocity  power  spectra  and  root-mean-square  (rms)  perttirbation  levels  were  taken  from 
the  experiment  of  Spencer  and  Jones^.  Because  phase  information  is  not  included  in  the  power 
spectra,  a  random  phase  relationship  for  the  modes  comprising  the  spectra  was  assumed.  The 
velocity  components  are  found  by  performing  a  Fourier  transform  of  the  complex  Fourier  coeffi¬ 
cients,  a,{0,  defined  by: 


0,(0  = 


(cosY  +  isiny) 


(5) 


where /is  the  frequency,  is  the  normalized  spectrum  function  measured  experimentally,  m7  is 

the  rms  level  of  the  velocity  component  also  measured  experimentally,  T  is  the  time  interval  to 
generate  boundary  conditions  (typically  four  times  the  total  period  of  simulation  desired),  and  y 
is  the  randomly  generated  phase  angle.  The  complete  details  of  the  derivation  of  time-dependent 
inlet  boundary  conditions  based  on  a  experimentally  measured  spectra  are  given  in  Wilson^. 
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The  convective  boundary  condition  is  applied  at  the  outflow  boundary  in  the  current  simula¬ 
tions  to  allow  the  developing  structures  to  flow  smoothly  out  of  the  computational  domain.  Writ¬ 
ten  in  terms  of  the  general  variable,  ^ ,  the  convective  boundary  condition  is  given  by. 


+  =  0 
3/  '3x, 


(6) 


No-stress  boundary  conditions  are  applied  at  the  free-streams  located  at  the  top  and  bottom  of 
the  computational  domain.  In  addition,  the  ffee-stream  boundary  is  placed  at  a  transverse  distance 
such  that  there  are  no  interactions  between  the  developing  structures  and  the  boundaries,  which 
would  violate  the  no-stress  condition.  The  no-stress  boundary  conditions  are  approximated  as  p  = 
constant  along  the  upper  and  lower  boundaries  and  du/dx^  =  0. 


The  cell-centered,  finite-volume  approach  is  used  to  integrate  the  continuous  governing  equa¬ 
tions,  given  by  equations  (1)  and  (2),  over  a  typical  nine-point  control  volume  which  yields  a  dis¬ 
cretized  system  of  algebraic  equations  to  be  solved  simultaneously.  Simulations  are  performed 
with  a  third-order  upwind  biased  approximation  of  convection  terms  following  Leonard^^.  The 
diffusive  flux  at  the  cell-face  is  approximated  using  fourth-order  central  differences.  Thus,  the  for¬ 
mal  spatial  accuracy  of  the  numerical  scheme  is  third-order.  A  third-order  fully  implicit  time  dis¬ 
cretization  is  applied  to  the  temporal  derivative.  Hence,  iteration  is  required  at  each  time  step. 
Details  of  the  numerical  scheme  are  presented  in  Wilson^.  A  typical  simulation  with  a  grid  of 
320  X  100  over  a  domain  size  of  L  x  //  =  20  jc  8,  requires  approximately  3  1/2  hours  of  CPU  time 
on  a  single  processor  of  the  Cray  C-90  supercomputer. 
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IIL  ACCURACY  OF  NUMERICAL  FORMULATION 

The  accuracy  of  the  numerical  method  is  checked  by  comparison  of  the  growth  rates  of  single 
frequency  disturbances  with  those  predicted  by  inviscid  theory.  Minimum  grid  points  for  accurate 
resolution  of  each  wave  can  thus  be  determined.  Given  a  base  velocity  profile,  linear  stability  the¬ 
ory  can  predict  the  initial  growth  rates  of  instability  waves.  At  high  Reynolds  numbers  viscous 
effects  are  negligible  so  that  inviscid  theory  can  be  used.  Monkewitz  and  Huerre^^  have  computed 
growth  rates  for  a  range  of  ^  (the  non-dimensional  angular  frequency,  =2nfQyUf.)  using  linear 
inviscid  stability  theory.  They  found  that  modes  with  p  <  0.5  had  exponentially  increasing  growth 
rates  whereas  those  with  p  >  0.5  experienced  exponentially  decreasing  growth  rates.  The  most 
amplified  mode  had  p  =  0.22  and  this  mode  was  referred  to  as  the  fundamental  mode.  Miksad^^ 
showed  that  linear  theory  is  appropriate  over  the  first  1  1/2  wavelengths.  Present  computations  of 
the  growth  rates  are  compared  to  the  results  of  Monkewitz  and  Huerre^^  in  Fig.  2  for  the  velocity 
ratio,  X  =  0.5 .  For  computations  in  Section  IV,  the  velocity  ratio,  X  =  0.538 ,  is  utilized.  Computa¬ 
tions  are  performed  at  Re  =  10^  for  two  perturbation  frequencies,  the  fundamental  (p  =  0.22)  and 
its  subharmonic  ;p  =  0.11 ),  with  the  initial  forcing  strength  for  both  modes  being,  = 

4.9  X  10"^.  Three  streamwise  grid  resolutions  are  utilized.  The  transverse  resolution  is  chosen  to 
resolve  the  eigenfunctions  of  the  two  modes  and  is  constant  for  aU  grids  used  in  this  section.  The 
information  on  the  grid  distribution  and  the  domain  size  is  summarized  in  Table  1.  Also  included 
is  the  numbo*  of  streamwise  points  per  wavelength  ippw).  For  the  fundamental  mode,  the  coarse 
grid  has  8  ppw  and  is  seen  to  be  inadequate  in  predicting  the  growth  rate.  The  intermediate  and 
fine  grids  have  16  ppw  and  32  ppw  respectively  and  are  able  to  predict  the  growth  rate  found  from 
inviscid  theory  to  within  roughly  2%.  For  the  subharmonic  mode,  the  coarse  (16  ppw),  intermedi¬ 
ate  (32  ppw),  and  fine  grids  (64  ppw)  predict  growth  rates  which  are  widiin  roughly  4%  of  that 
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predicted  by  linear  theory.  Therefore,  with  the  present  numerical  method,  at  least  16  ppw  are 
required  to  accurately  resolve  each  mode.  A  comparison  of  the  computed  growth  rates  and  those 
from  linear  theory  is  given  in  Table  2.  The  simulation  at  Re  =  10^  is  essentially  the  same  as  those 

at  Re  =  10^,  whereas  at  Re  =  10^,  the  predicted  growth  rales  (not  shown)  are  about  20%  lower. 

Table  1:  Domain  sizes  for  spatial  growth  rate  computations. 


Mode 

Coarse, 
ni  X  nj,  ppw 

Intermediate, 
ni  X  nj,  ppw 

Fine, 

ni  X  nj,  ppw 

LxH 

first  subharmonic 

48x80, 16 

96x80, 32 

192x80, 64 

5.84x4 

fundamental 

24x80,8 

48x80, 16 

96x80, 32 

2.92x4 

Table  2:  Comparison  of  growth  rates  from  simulation  to  linear  theory,  Re  =  10'* 


Mode 

Coarse 

Intermediate 

Fine 

Linear  theory 

first  subharmonic 

0.071 

0.071 

0.071 

0.074 

fundamental 

0.085 

0.096 

0.096 

0.098 

IV.  RESULTS  AND  DISCUSSION 

The  simulations  forced  with  discrete  modes  can  be  viewed  as  unsteady,  laminar  calculations 
of  instability  modes,  while  those  which  use  a  broad-band  spectrum  as  forcing  functions,  can  be 
viewed  as  a  2-D  approximation  of  the  3-D  turbulent  mixing  layer  observed  experimentally^^.  Fig. 
3  compares  the  spectrum  from  a  3-D  experiment  to  those  predicted  in  the  present  simulations 
using  (0  the  fundamental  mode  and  (i'O  a  broad-band  of  modes.  The  3-D  experimental  spectrum 
shows  a  much  more  gradual  decay  of  the  amount  of  energy  in  the  high  frequency  range  than 
obtained  in  our  simulation  with  broad  mode  forcing.  In  the  latter,  the  energy  is  concentrated  in  a 
broad  band  of  modes,  the  highest  frequency  modes  having  been  removed  by  dissipation,  but  the 
2-D  approximation  prevents  replenishment  via  the  vortex  stretching  mechanism.  Results  of  simu- 
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lations  forced  with  a  single  discrete  mode  show  a  concentration  of  energy  at  that  associated  fre¬ 
quency,  as  expected. 

The  resulting  vortex  dynamics  for  mixing  layer  simulations  using  the  two  sets  of  perturbations 
at  the  inlet  boundary  are  presented  in  this  section.  The  effect  of  the  vortex  dynamics  on  time-aver¬ 
aged  properties  are  discussed  in  the  next  section.  Effects  of  inlet  forcing  and  Reynolds  numbers 
are  explored. 

When  no  forcing  is  applied  at  the  inlet  boundary  (i.e.  steady  inlet  boundary  conditions),  time 
traces  of  velocity  and  pressure  confirm  that  little  unsteady  motion  occurs  after  initial  transients 
have  been  convected  out  of  the  computational  domain.  The  maximum  normalized  fluctuating 
velocity  within  the  computational  domain  was  found  to  be  =  0.001,  where  is 

the  rms  of  the  streamwise  fluctuating  velocity  and  Am  =  f///  -  Ui  If  some  form  of  forcing  is  used 
at  the  inlet  boundary  for  the  period,  0  <  t  <  x,  after  which  the  forcing  is  suspended  for,  r  >  x,  the 
unsteady  motion  due  to  the  initial  forcing  is  convected  out  of  the  computational  domain  and 
steady  flow  results.  This  finding  is  in  agreement  with  the  spatial  simulations  of  Lowery  and 
Reynolds^'  in  that  unsteady  motion  would  develop  from  steady  inlet  conditions  only  for 
extremely  long  computational  domains  and/or  very  large  Reynolds  numbers.  However,  such 
steady  inlet  conditions  are  impossible  to  achieve  in  physical  experiments.  Some  form  of  forcing  is 
needed  to  initiate  the  instabilities  of  the  mixing  layer  at  moderate  Reynolds  numbers  and  realistic 
computational  lengths. 
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Vortex  Dynamics 

Forcing  with  Modes  from  Inviscid  Stability  Theory 

In  this  section,  inlet  boundary  perturbations  are  derived  from  the  inviscid  stability  analysis  of 
Monkewitz  and  Huerre^^.  Some  aspects  of  discrete  mode  simulations  were  discussed  in  the  previ¬ 
ous  section  where  single  mode  simulations  were  utilized  to  test  the  grid  resolution  requirements 
of  the  numerical  formulation.  A  more  detailed  picture  of  these  simulations  is  given  in  this  section 
as  well  as  the  results  from  two  and  three  mode  simulations. 

The  time,  t,  is  non-dimensionalized  by  the  eddy  flow-through  time,  t,  =  L/Uc ,  which  is  the 
time  required  for  a  structure  to  travel  from  the  inflow  to  the  outflow  at  the  convective  velocity. 
The  vorticity  contours  at  t  =  2.708  seconds,  for  a  simulation  forced  by  the  fundamental  mode,  are 
shown  in  Fig.  4  for  a  computational  domain  of  L  x  //  =  12  x  4  with  ni  x  nj  =  200  x  100  grid  points. 
This  simulation  shows  that  the  initially  uniform  vortex  sheet  at  /  =  0,  “rolls  up”  into  intense 
“lumps”  of  spanwise  vorticity  separated  by  a  distance  of  the  fundamental  wavelength  (A,  =  0.985). 
This  process  is  the  so-called  primary  or  Kelvin-Helmholtz  instability.  After  a  distance  of  roughly 
2  1/2  wavelengths  (x  =  2.5)  from  the  inflow,  the  roll-up  of  the  shear  layer  saturates  and  experi¬ 
ences  a  decay  for  the  remaining  spatial  extent. 

Many  experimental  and  numerical  studies  have  shown  that  in  the  absence  of  subharmonic  fre¬ 
quencies,  vortex  pairing  does  not  occur.  To  test  the  effects  of  dual-mode  forcing,  simulations  were 
carried  out  with  forcing  consisting  of  the  fundamental  and  first  subharmonic  modes  at  equal  initial 
amplitudes  such  that  the  rms  of  streamwise  velocity  at  the  inlet  was  {m,„/Am)^^^  =  0.06.  This 
simulation  was  carried  out  using  the  same  computational  grid  and  duration  as  the  fundamental 
mode  simulation  previously  discussed.  Vorticity  contours  for  the  two-mode  sim  tion  at  /  = 
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2.708  seconds  are  shown  in  Fig.  5.  After  the  initial  start-up  period,  regular  periodic  motion  occurs 
for  t  >  1.625  seconds.  The  first  two  fundamental  wavelengths  (0  <  or  <  2)  are  similar  to  the  single 
mode  simulation,  while  downstream  the  flow  is  drastically  different.  At  x  =  4,  the  subharmonic 
mode  is  of  sufficient  strength  to  vertically  displace  the  lumps  of  vorticity  separated  by  fundamen¬ 
tal  wavelengths.  This  causes  the  two  neighboring  vortices  to  rotate  about  their  common  renter 
and  coalesce  forming  a  new  larger  structure  which  results  in  a  rapid  growth  of  mixing  la>  th 
at  the  location  of  vortex  pairing.  Downstream  of  the  pairing  location,  the  new  structure  rotates 
clockwise  through  a  small  angle  and  is  convected  out  of  the  domain  with  the  spacing  between  the 
new  structures  being  one  subharmonic  wavelength  (A  =  1.97). 

A  three-mode  simulation  was  carried  out  using  the  fundamental,  first,  and  second  subhar¬ 
monic  modes  from  inviscid  theory  as  inlet  forcing  functions  such  that  =  0.05.  The 

computational  length  and  height  of  the  domain  were  increased  to  L  x  //  =  20  x  8,  with  ni  x  nj  =  320 
X  100  grid  points.  An  expanding  grid  was  used  in  the  transverse  direction  to  enable  the  rapid 
changes  in  the  fundamental  mode  eigenfunctions  near  the  center  of  the  mixing  region  to  be 
resolved  efficiently.  One  period  of  the  second  subharmonic  mode  is  shown  in  Fig.  6,  in  which 
only  a  section  of  the  computational  height  (2  <  y  <  6)  is  displayed.  Note,  that  the  structures  move 
with  the  convective  speed  of  the  flow,  Uc  =  0.65,  therefore  the  structures  advance  approximately 
1.3  units  in  the  downstream  direction  in  the  plotting  interval  (0.065  seconds)  shown  in  Fig.  6.  The 
addition  of  the  second  subharmonic  mode  allows  a  second  vortex  pairing  to  occur  at  x  =  10.  The 
time  trace  and  power  spectra  results  (not  shown)  reveal  the  transfer  of  energy  from  the  fundamen¬ 
tal  to  the  first  subharmonic  mode  (the  first  vortex  pairing)  and  then  from  the  first  to  the  second 
subharmonic  mode  (the  second  vortex  pairing)  as  the  structures  convect  from  inflow  to  outflow. 
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Forcing  with  a  Broad  Spectrum  of  Modes  from  Experimental  Data 


The  second  class  of  perturbations  used  as  forcing  functions  is  that  of  a  broad  spectrum  of 
modes  derived  from  experimentally  measured  power  spectra.  This  class  of  inlet  boundary  condi¬ 
tions  is  used  to  simulate  the  broad  range  of  scales  found  in  fully  three-dimensional  turbulent 
flows.  Of  particular  interest  is  the  pairing  mechanism  and  shear  layer  growth  rates  compared  to 
one-,  two-,  and  three-mode  simulations. 

Many  con^utational  and  experimental  mixing  layer  studies  have  shown  that  small  scale 
motions  have  their  origins  in  the  braid  region  of  the  two-dimensional  structures  where  streamwise 
vorticiQ^  first  appears.  Obviously,  2-D  simulations  are  not  capable  of  capturing  the  energy  cascade 
to  small  scales  and  can  only  simulate  the  generation,  growth,  and  interaction  of  the  2-D  spanwise 
vortex  structures.  The  current  simulations  confirm  that  for  2-D  simulations,  small  scale  motions 
are  not  continuously  generated.  Moreover,  if  small  scale  motions  are  injected  at  the  inlet  bound¬ 
ary,  they  experience  rapid  decay  in  agreement  with  linear  theory  as  discussed  previously. 

It  has  been  shown,  using  discrete  mode  simulations,  that  periodic  forcing  with  no  explicit 
phase  shift  between  the  modes  produces  complete  vortex  pairing  at  fixed  spatial  locations.  The 
non-periodic  inlet  boundary  forcing  used  with  the  broad  mode  simulations  does  not  produce  vor¬ 
tex  pairing  at  fixed  spatial  locations  and  pairing  occurs  somewhat  randomly  over  a  region.  Vortex 
tearing,  whereby  a  weaker  vortex  is  pulled  apart  by  its  stronger  neighboring  vortices,  can  be  seen 
(not  shown)  in  selected  vorticity  contours  of  broad  mode  simulations,  due  mainly  to  the  presence 
of  finite  phase  shifts  between  the  various  modes  in  the  inlet  forcing  functions.  The  effect  of  the 
strength  of  the  forcing  is  compared  in  Hg,  7a  -  7c,  which  shows  vorticity  contours  at  r  =  3.413 
seconds  and  Re  =  10^  for  forcing  levels  of  =  0.0015,  0.015, 0.15,  respectively.  The 


simulations  were  carried  out  with  the  same  computational  domain  and  grid  as  was  used  in  the 
three  mode  simulation.  With  the  full  strength  forcing  amplitude  reduced  by  a  factor  of  100,  Fig. 
7a,  shear  layer  roll-up  occurs  at  jf  =  10,  but  there  is  no  evidence  of  pairing  within  the  computa¬ 
tional  domain,  0  <x  <20.  \Wth  the  forcing  amplitude  reduced  by  a  factor  of  10,  Fig.  7b,  shear 
layer  roll-up  occurs  at  around  x  =  6,  with  the  first  pairing  around  x  =  14.  The  strongly  forced  case. 
Fig.  7c,  shows  that  shear  layer  roll-up  occurs  immediately  downstream  of  the  inlet  boundary  and 
the  first  pairing  occurs  at  around  x  =  4.5.  These  trends  are  repeated  in  Fig.  8,  for  simulations  at  Re 
=  10^.  The  effect  of  forcing  strength  on  time-averaged  turbulent  statistics  is  discussed  in  a  later 
section. 

Clearly,  the  effect  of  inlet  forcing  is  to  promote  rapid  growth  of  instabilities  to  amplitudes 
which  are  sufficiently  high  for  shear  layer  roll-ups  to  occur,  and  in  the  presence  of  sub-harmonics 
vortex-pairings  would  take  place  subsequently.  The  far-field  vortex  dynamics  are  mostly  indepen¬ 
dent  of  the  inlet  forcing.  The  fluctuating  intensity  of  0.15%  corresponding  to  the  forcing  levels  in 
Figures  7a  and  8a  is  of  the  order  of  the  background  turbulence  in  many  wind  tunnels.  In  fact,  a 
simulation  with  no  forcing,  i.e.,  completely  steady  inlet  conditions  deviates  more  from  practical 
experimental  conditions  than  one  with  small  forcing.  One  may  speculate  that  given  sufficient 
length,  in  the  former  case,  minute  instabilities  in  the  form  of  round-off  errors  in  the  computer 
would  ultimately  grow  sufficiently  to  produce  the  same  vortex  dynamics.  Needless  to  say,  there 
has  been  no  simulation  of  sufficient  length  to  test  this  hypothesis. 

To  confirm  that  the  observed  vortex  dynamics  would  continue  downstream  a  simulation  was 
carried  out  for  an  expanded  domain  size  ofLxf/  =  40xl6  with  ni  x  nj  =  640  x  130  grid  points. 
Figures  7d  and  8d  show  the  resulting  vorticity  contours  at  f  =  1.625  seconds  for  Re  =  10^  and  10^, 
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respectively,  using  the  broad  mode  forcing  with  the  amplitude  (Urm/^Kax  ~  0  015.  The  results 
at  this  time  step,  as  well  as  others  not  shown,  show  that  the  mixing  layer  continues  to  grow 
through  vortex  pairing  for  the  additional  length,  20  <  x  <  40.  The  time  trace  results  at  various 
points  along  the  length  of  the  layer  centerline  show  that  the  dominant  frequency  shifts  from  /  = 
0.8  at  jc  =  3  to  /  =  0.08  at  jc  =  36.  It  is  apparent  that  on  the  average,  approximately  four  pairings 
have  taken  place  because  the  dominant  frequency  has  been  halved  four  times. 

Forcing  with  Mixed  Modes 

To  approximate  experiments  in  which  the  shear  layer  is  forced  with  a  few  discrete  modes,  two 
additional  simulations  at  Re  =10^  are  discussed  in  this  section.  The  first  simulation,  which  will  be 
referred  to  as  case  b,  uses  the  fundamental  mode  with  strength,  =  0.015,  superim¬ 

posed  on  the  broad  mode  forcing  function  used  in  the  weakly  forced  simulation  (Fig  8a,  referred 
to  as  case  a),  which  had  a  strength,  =  0.(X)15.  The  second  simulation,  refeired  to  as 

case  c,  uses  the  fundamental  and  the  first  subharmonic  modes  of  equal  strength,  = 

0.(X)15,  along  with  the  weak  broad  modes  used  in  case  b.  Vorticity  contours  at  r  =  3.413,  are  com¬ 
pared  in  Rgures  9a  -  9c  for  cases  a  -  c,  respectively.  The  addition  of  the  fundamental  mode  which 
is  one  order  of  magnitude  larger  than  the  broad  mode  forcing  function  (Fig  9b)  results  in  shear 
layer  roll-up  at  an  earlier  streamwise  location  ofx  =  2.  However,  the  location  of  the  initial  vortex 
pairing  is  shifted  downstream  from  x  =  8.5  (Fig  9a)  to  x  =  9.5  (Fig  9b).  The  addition  of  the  funda¬ 
mental  and  first  subharmonic  modes  (Fig  9c)  results  in  shear  layer  roll-up  and  a  large  shift  in  the 
location  of  the  first  vortex  pairing,  from  x  =  8.5  (Fig  9a)  to  x  =  4  (Fig  9c),  but  with  no  vortex  inter¬ 
action  in  the  region,  5  <  x  <  13.  Hme  trace  results  (not  shown)  indicate  that  the  second  subhar¬ 
monic  mode  which  is  present  in  the  weak  broad  mode  forcing  increases  gradually  in  this  region, 
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but  does  not  grow  to  sufficient  strength  to  produce  a  second  pairing  within  the  computational 
domain,  unlike  in  the  simulations  of  cases  a  and  b.  It  appears  that  the  addition  of  the  first  subhar- 
monic  mode  in  case  c  suppresses  the  growth  of  the  second  subharmonic  mode*^  thereby  resulting 
in  a  large  region  with  no  vortex  interaction  and  stagnation  of  mixing  layer  growth.  The  momen¬ 
tum  thickness  of  the  shear  layer  is  defined  by  the  following  integral; 

H 

0  =  («H  -«)(«-  u^)  dy  (7) 

0 

Hme-averaged  momentum  thicknesses  for  the  three  cases  are  shown  in  Fig  9d.  It  can  be  seen  that 
the  addition  of  the  fundamental  mode  in  case  b,  results  in  a  slight  increase  in  layer  width  in  the 
region,  2  <  jc  <  4,  followed  by  a  slight  decrease  in  the  region,  5  <  x  <  10.  The  addition  of  the  ftin- 
damental  and  first  subharmonic  modes  in  case  c,  results  in  a  significant  increase  in  shear  layer 
width  in  the  region,  1  <  x  <  6,  followed  essentially  by  no  growth  in  the  region,  6  <  x  <  14.  Shear 
layer  growth  is  resumed  for  x  >  14  with  the  growth  of  the  second  subharmonic  mode.  From  the 
previous  discussion,  the  role  of  the  low  amplitude  broad  mode  forcing  is  to  provide  a  supply  of 
suUiarmonics  to  produce  vortex  pairing  after  the  discretely  forced  modes  have  saturated.  This  is 
apparent  by  comparing  the  simulation  shown  in  Fig  4  (fundamental  mode  only)  to  that  shown  in 
Fig  9b  (fundamental  mode  and  weak  broad  mode  forcing). 

Three-dimensional,  temporal  simulations  utilizing  a  discrete  mode  superimposed  on  a  broad 
spectrum  of  modes  were  performed  by  Riley  and  Metcalfe  Their  results  show  that  the  addition 
of  the  fundamental  mode  to  strong  background  turbulence  (of  18%  intensity),  which  is  similar  in 
spectrum  to  the  In-oad  band  forcing  used  in  this  paper,  resulted  in  an  increased  initial  growth  of  the 
shear  layer,  followed  by  a  suppressed  growth  rate  in  comparison  to  results  using  only  background 
turbulence.  Those  trends  are  in  agreement  with  the  results  of  case  a  and  b.  The  experimental  study 
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of  fwced  mixing  layers  by  Oster  and  Wygananski^^  utilized  a  single  frequency  disturbance  to 
force  the  shear  layer.  Their  results  showed  that  the  unforced  mixing  layer  (with  only  a  small  level 
of  background  turbulence  of  roughly  0.2%  intensity)  exhibited  constant  linear  growth,  while  the 
addition  of  the  discrete  mode  showed  increased  growth  rates  in  the  near  field,  followed  by  regions 
of  zero  growth  in  which  it  was  shown  that  no  vortex  interaction  took  place,  followed  by  resumed 
mixing  layer  growth.  The  present  results  are  in  complete  agreement  with  those  experimental 
observations. 

Effect  of  Reynolds  Number 

The  effect  of  Reynolds  number  for  a  single  mode  simulation  has  already  been  discussed.  It 
was  shown  that  the  medium  to  high  Reynolds  number  simulations  {Re  >  10^)  yield  close  agree¬ 
ment  with  spatial  growth  rates  predicted  by  linear  inviscid  theory.  For  low  Reynolds  number  sim¬ 
ulations,  Re^  =  100,  the  spatial  growth  rates  are  mildly  damped  compared  to  linear  inviscid 
theory  and  there  is  a  mild  smearing  of  the  rolled-up  structures.  Davis  and  Moore^^  report  similar 
trends  in  their  one-,  two-,  and  three-mode  simulations,  with  reduced  Reynolds  number  simula¬ 
tions  affecting  only  the  sharpness  of  the  vortical  structures.  The  authors  report  that  the  vortex 
dynamics  are  not  changed  by  reducing  the  Reynolds  number  from  10^  to  25.  These  results  support 
the  claim  that  for  two-dimensional  simulations  involving  only  a  few  discrete  modes  from  linear 
stability  theory  as  forcing  functions,  the  Reynolds  number  is  of  minor  importance  in  computing 
the  vortex  dynainics  of  the  shear  layer. 

The  broad  spectrum  forcing  function  contains  many  high  frequency  modes  which  are  more 
susceptible  to  viscous  forces  than  their  low  frequency  counterparts.  Simulations  using  the  broad 
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spectrum  forcing  function  were  carried  out  aX  Re  =  lO'*  to  determine  if  the  decreased  damping  of 
high  frequency  modes  would  result  in  a  change  in  the  pairing  mechanism  discussed  previously. 

An  inspection  of  the  vorticity  contours  in  Fig.  8  for  the  /?e  =  10^  simulations  reveals  that  the 
shear  layer  rolls  up  and  pairs  much  earlier  in  comparison  to  the  Re  -  10^  simulations.  Also,  the 
vortices  tend  to  maintain  their  individual  identities  after  pairing  while  the  increased  viscous  forces 
of  the  Re  =  10  simulations,  tend  to  smear  out  the  individual  vortices  and  create  a  more  uniform 
core  for  the  paired  structures. 

The  effect  of  the  strength  of  forcing  is  compared  in  Fig.  8  for  the  Re  -  10^  simulations.  The 
forcing  levels  in  Figures  8a  -  c  are  =  0.0015,  0.015,  and  0.15,  respectively,  with 

identical  inlet  boundary  conditions  and  grid  used  for  the  Re  =  10^  simulation  shown  in  Fig.  7. 
\A^th  the  smallest  forcing  amplitude.  Fig.  8a,  the  shear  layer  rolls  up  at  around  x  =  5  and  the  first 
pairing  occurs  at  around  x  =  9.  Recall  that  the  same  forcing  level  produced  no  vortex  pairing  for 
the  Re  =  10  simulation.  With  the  intermediate  forcing  amplitude.  Fig.  8b,  shear  layer  roll-up 
occurs  at  about  x  =  3,  with  the  first  pairing  at  around  x  =  6.  The  full  strength  forcing  results,  shown 
in  Fig.  8c  show  that  the  layer  rolls  up  very  close  to  the  inflow  boundary  and  the  first  pairing 
occurs  around  x  =  3.  This  confirms  the  conclusions  reached  using  single  mode  simulations,  in  that 
decreasing  the  Reynolds  number  has  a  damping  effect  on  the  growth  rate  of  disturbances.  The 
simulations  at  /?e  =  10^  suggest  that  the  streamwise  location  of  the  initial  vortex  pairing  scales  lin¬ 
early  with  the  log  of  the  inlet  forcing  amplitude,  which  is  not  surprising  since  linear  theory  pre¬ 
dicts  exponential  growth  rates  for  these  instability  waves. 


20 


Time-Averaged  Properties 


This  section  discusses  the  time-averaged  properties  of  the  mixing  layer.  In  particular,  the 
effect  of  the  inlet  boundary  conditions  and  coherent  structures  on  mixing  layer  growth,  time-aver¬ 
aged  velocity  profiles,  and  turbulence  statistics  are  investigated. 

Discrete- Mode  Forcing 

Three  normalized  shear  layer  widths  are  calculated  with  the  first  based  on  the  transverse  dis¬ 
tance  between  locations  where  the  property  ^  =  0.9  and  ^  =0.1,  with  (|>  =  («  -  /  (««  -  «£,)  •  The 

second  normalized  width  is  the  momentum  thickness  which  was  defined  previously  The  third 
width  is  the  normalized  vorticity  thickness,  5^,  which  was  defined  previously  in  Section  II.  The 
vorticity  thickness  shows  more  drastic  changes  in  shear  layer  width  because  it  is  a  differential 
quantity,  while  the  momentum  thickness  is  an  integral  quantity  and  is  thus  smoother. 

Figures  10  and  1 1  show  time-averaged  results  for  the  three  mode  simulation  averaged  over 
three  periods  of  the  second  subharmonic  mode.  The  shear  layer  widths  displayed  in  Fig.  10a 
clearly  show  the  effect  of  the  two  spatially  fixed  vortex-pairing  events  which  occur  roughly  at  x  = 
4  and  x  =  10.  The  shear  layer  width  experiences  a  drastic  reduction  after  the  second  vortex  pairing 
around  the  region  x  =  13.5,  after  which  the  layer  growth  is  arrested.  The  cause  of  the  decrease  in 
layer  width  can  be  determined  by  examining  the  vorticity  contours  shown  in  Fig.  6.  The  second 
vortex  pairing  occurs  atx  =  10  where  the  two  vortices  are  aligned  vertically  and  reach  maximum 
transverse  separation.  Downstream  of  this  point,  the  elongated  structure  rotates  clockwise  and 
achieves  a  horizontal  orientation  at  roughly  x  =  13.5  (seen  at  f  =  1.755  seconds,  Fig.  6c).  Still  fur¬ 
ther  downstream,  the  structure  thickens  in  the  transverse  direction  slightly  as  it  convects  through 
the  outflow  boundary.  The  horizontal  orientation  of  the  elongated  structure  at  x  =  13.5  results  in  a 
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decrease  in  mixing  layer  width  along  with  negative  Reynolds  shear  stress.  The  resulting  negative 
Reynolds  shear  stress  can  be  deduced  from  the  orientation  of  the  vortex  pairing  as  noted  by  Ho 
and  Huerre^^.  Before  the  two  vortices  pair  and  reach  their  vertical  orientation  (Fig.  6c,  jc  =  9.5), 
the  correlation  of  vertical  and  horizontal  velocity  fluctuations  produces  positive  Reynolds 
stresses,  while  after  they  reach  their  vertical  orientation  (Fig.  6d,  x  =  11),  negative  Reynolds  shear 
stress  results. 

Fig.  10b  displays  the  time-averaged  u  component  profiles  which  show  self-similarity  except  at 
cross  sections  where  pairings  occur  (x  =  4  and  x  =  10).  The  v  component  profiles,  after  the  second 
pairing  atx  =  12.4,  shown  in  Fig.  10c,  reveal  negative  mean  transverse  velocity  at  the  bottom  of 
the  computational  domain  which  results  in  transpoit  of  fluid  away  from  the  mixing  layer.  Fig.  1  la 
shows  characteristic  double  “hump”  profiles  for  cross  sections  after  the  second  pairing  (x  > 
12)  and  Fig.  11c  shows  relatively  large  negative  Reynolds  shear  stress  for  the  region  of  decreasing 
layer  width  (x  =  12.4)  as  explained  above.  The  turbulence  statistics  are  compared  to  experimental 
data  from  the  forced  shear  layer  experiments  of  Oster  and  Wygnanski^^.  The  double  hump 
profile  is  also  apparent  in  this  experiment,  although  it  is  not  as  pronounced  as  in  the  current  2-D 
simulations. 


Broad-Spectrum-Mode  Forcing 

Recall  that  forcing  with  periodic  inlet  boundary  conditions  from  inviscid  theory  produces  vor¬ 
tex  pairing  at  fixed  spatial  locations,  while  quasi-random  broad  spectrum  forcing  produces  vortex 
pairings  which  are  distributed  over  a  range  of  spatial  locations.  As  a  result,  the  time-averaged  vor- 
ticity,  pressure,  and  shear  layer  widths  vary  more  smoothly  than  the  time-averaged  properties 
from  the  three  mode  simulation.  The  broad  mode  simulations  at  Re  =  10^  and  Re  =  lO'^,  Figures 


22 


7d  and  8d,  respectively,  are  used  to  calculate  spread  rates,  mean  velocity  profiles,  and  turbulence 
statistics  over  the  interval,  1.056  <  t  <  4.063  seconds  which  represents  roughly  3  eddy  flow¬ 
through  times.  From  analysis  of  time  trace  results  of  the  Re  =  10^  simulation,  roughly  120  struc¬ 
tures  have  been  generated  just  downstream  of  the  inflow  and  about  20  structures  have  passed 
through  the  outflow  in  this  time  interval.  The  time-averaged  results  for  these  simulations  are  pre¬ 
sented  in  Fig.  12  and  13.  In  comparison  to  the  3'mode  simulation  (Fig.  10a),  the  continuous  sup¬ 
ply  of  subharmonics  in  this  case  results  in  continued  linear  growth  of  the  shear  layer.  Mean 
streamwise  velocity  profiles  show  self-similarity  throughout  the  domain,  while  naean  transverse 
velocity  profiles  show  that  on  the  average  fluid  is  entrained  into  the  mixing  region.  Profiles  of 
streamwise  fluctuating  velocities  (Fig.  13a)  show  that  the  double  hump  profiles  found  in  the  3- 
mode  simulations  are  absent  The  results  are  in  reasonable  agreement  with  experimental  observa¬ 
tions  for  unforced  shear  layers.  Reynolds  shear  stress  profiles  (Fig.  13c)  are  positive  and  similar  in 
magnitude  to  those  measured  in  the  vortex  interaction  region  of  Oster  and  Wygnanski^^.  The 
magnitude  of  the  transverse  fluctuating  velocity  is  overpredicted  by  approximately  40%  in  the 
current  2-D  simulations.  Lowery  and  Reynolds^^  made  similar  observations  in  their  2-D  simula¬ 
tions  while  their  3-D  simulations  yielded  the  experimentally  measured  magnitude.  The  discrep¬ 
ancy  is  clearly  due  to  the  absence  of  energy  transfer  to  small  scales  due  to  vortex  stretching.  The 
sante  general  trends  are  observed  at  both  Reynolds  numbers,  except  that  the  profiles  vary  more 
sharply  at  the  higher  Reynolds  number.  Spread  rates  computed  from  the  discrete  and  broad  spec¬ 
trum  simulations  are  presented  and  compared  to  experimental  values  in  Table  3  in  the  following 
normalized  form: 


-L^  =  2X^ 

Audi  dx 


(8) 
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It  can  be  seen  that  the  spread  rate  in  the  three-mode  simulation  is  greatly  enhanced  in  the  region 
prior  to  the  second  vortex  pairing  event  (150  <  x/8  <  240).  Also,  the  broad  spectrum  simulations 
at  Re  =  10  and  10  yield  nearly  identical  spread  rates  in  the  far-field  which  are  similar  to  the  far- 
field  spread  rates  measured  experimentally. 


Table  3:  Comparison  of  computed  spread  rates  with  experimental  values. 


Study 

Re 

Region 

Spread  Rate, 

Current,  3-mode, 

IHHCSHI 

0<x/e  <120 

(Fig.  10a) 

130<Jc/e  <  140 

150<x/e  <240 

350  <  x/8  <  600 

Current,  broad 

10^ 

15O<x/0  <  1000 

U.022 

spectrum,  (Fig.  12a) 

Current,  broad 

10" 

15O<x/0  <  1000 

0.021 

spectrum,  (Fig.  12d) 

Wygnanski  & 

10^ 

- 

0.022 

Fieldler^^ 

Spencer  and  Jones^^ 

10^ 

- 

0.016 

Oster  and 

10^ 

0<x/d  <400 

Wygnanski^^ 

400<jc/e  <600 

600<j:/e  <  1500 

The  effect  of  the  forcing  strength  on  the  streamwise  development  of  the  maximum  cross-sec¬ 
tional  value  of  and  is  presented  in  Figures  14  and  15.  The  results  at  Re  =  10^  are  not 
completely  independent  of  forcing  level,  while  the  Re  =  lO'^  far-field  results  (jc>  10)  show  only 
minor  sensitivity  to  initial  forcing  level.  One  comes  to  the  conclusion  that  the  larger  growth  rate 
with  increased  Reynolds  number  results  in  asymptotic  values  at  earlier  streamwise  locations.  The 
computational  length,  L  =  20,  would  need  to  be  increased  to  achieve  similar  forcing  insensitivity 
in  the  /?^  =  10^  simulation.  Fig.  15  directly  compares  the  Re  =  10^  and  10^  simulations  for  the 
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intermediate  forcing  level  and  increased  length,  L  =  40.  With  the  length  doubled,  the  turbulent  sta 
tistics  from  the  Re  =  10^  simulation  slowly  approach  those  of  the  Re  =  10"^  simulation  in  the  far- 
field.  The  Ur^  values  appear  to  approach  the  asymptotic  values  more  slowly  than  the  values. 

Estimation  of  Length  and  Time  Scales 

IVo-point  spatial  correlations  and  autocorrelations  were  computed  at  locations  along  the  cen¬ 
terline  of  the  mixing  layer  to  estimate  the  size  and  lifetime  of  the  coherent  structures.  The  turbu¬ 
lent  statistics  of  the  plane  mixing  layer  are  not  isotropic,  therefore  streamwise  and  transverse 
correlations  are  considered  separately.  It  has  been  shown  that  the  coherent  structures  move  down¬ 
stream  at  the  convective  velocity.  Davies  et  al.^®  used  this  information  and  showed  that  an  enve¬ 
lope  of  maximum  correlation  can  be  formed  by  computing  the  autocorrelation  in  a  moving 
reference  frame  which  travels  with  the  convective  velocity  along  the  mixing  layer  centerline.  The 
autoconelation  is  then  given  by: 


H  i  0  Uj  (Jq  +  yo>  t + T) 


(9) 


where  Xg,  yg  is  the  location  where  the  correlation  is  computed  (y^  is  located  transversely  on  the 
layer  centerline  with  (|)  =  0.5),  and  x  is  the  temporal  separation  normalized  by  Au  and  .  Thus, 
the  reference  frame  moves  a  distance  of  r  =  UgX  along  the  layer  centerline  in  the  time,  x .  The 
autocorrelation  was  used  to  estimate  the  lifetime  of  the  coherent  structures  at  various  streamwise 
locations. 


The  /?ii  and  R22  autocorrelations  at  4  streamwise  locations  are  shown  in  Fig.  16.  The  results 
show  that  the  structures  at  the  first  location  (Fig.  1 6a  and  16e)  have  short  lifetimes  and  become 
decorrelated  in  a  relatively  short  period  of  time  compared  to  later  spatial  locations.  Thus,  the  iden- 
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tides  are  lost  through  the  vortex-pairing  mechanisms,  whose  wavelengths  increase  downstream. 
The  intersection  of  a  parabola,  fit  to  the  data  at  t  =0,  and  the  horizontal  axis  gives  an  estimate  of 
the  micro  time  scale  or  lifetime  of  the  coherent  structures.  These  results  do  not  confirm  the  pres¬ 
ence  of  a  feedback  mechanism  suggested  by  Dimotakis  and  Brown^  which  would  ensure  that  the 
structures  are  correlated  all  the  way  upstream  to  the  inlet  plane.  The  integral  time  scales  for  the  Re 
=  10^  simulation  are  roughly  50%  larger  than  those  for  the  Re  =  10^  simulation.  The  larger  growth 
rate  of  disturbances  at  Re  =  10"^  results  in  a  faster  turnover  time  with  the  structures  becoming 
decorrelated  at  earlier  times. 

The  transverse  size  of  the  coherent  structures  was  estimated  using  two-point  spatial  correla¬ 
tions  which  requires  prescribing  the  spatial  separation,  r,  along  the  transverse  direction  with  zero 
time  separation.  The  two-point  correlation  is  given  by: 


(10) 


The  /?22  spatial  correlations  are  shown  in  Fig.  17.  As  expected  the  transverse  size  of  the  struc¬ 
tures  increase  from  inflow  to  outflow,  with  similar  trends  for  both  Reynolds  number  simulations. 


V.  CONCLUSION 

Numerical  simulations  of  the  spatially-developing  two-dimensional  mixing  layer  have  been 
performed  with  a  variety  of  unsteady  perturbations  of  the  inlet  boundary  conditions  used  as  forc¬ 
ing  functions.  The  forcing  functions  are  used  to  trigger  the  instabilities  which  occur  naturally  in 
mixing  layer  flows.  Moreover,  these  instabilities  can  be  exploited  to  alter  drastically  the  time- 
averaged  properties  of  the  flow. 
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Two  types  of  perturbations  were  considered;  (i)  one-,  two-,  and  three-mode  combinations  of 
eigenfunctions  corresponding  to  the  solution  of  the  inviscid  stability  equations  and;  (ii)  random 
inlet  boundary  conditions  derived  from  experimentally  measured  power  spectra.  The  simulations 
using  a  single  mode  from  inviscid  theory  show  that  over  the  first  wavelength  the  perturbations 
experience  exponential  growth  resulting  in  shear  layer  roll-up,  followed  by  saturation  and  even¬ 
tual  decay.  The  presence  of  subharmonic  modes  in  the  inlet  boundary  conditions  produced  vortex 
pairing,  with  the  number  of  such  pairing  events  equal  to  the  number  of  subharmonics  present  in 
the  inlet  boundary  conditions.  Saturation  occurred  once  the  supply  of  subharmonics  were 
exhausted  followed  by  decay.  Thus,  enhanced  mixing  due  to  discrete  mode  forcing  occurs  only 
over  limited  distances.  For  simulations  forced  with  a  broad  spectrum  of  modes  derived  from 
experiment,  vortex  pairing  occurred  over  a  region  but  not  at  fixed  locations  as  with  the  discrete¬ 
mode  simulations.  Three  vortex  mergings  at  one  spatial  location  could  be  observed  in  the  broad 
spectrum  simulations  as  weU  as  the  process  of  vortex  shredding. 

The  Reynolds  number  was  shown  to  greatly  affect  the  locations  of  the  initial  shear  layer  roU- 
up  in  the  broad  mode  simulations.  For  low  Reynolds  number  simulations,  the  additional  viscous 
forces  damp  the  development  of  shear  layer  roll-up  resulting  in  delayed  roll-up  and  vortex  pair¬ 
ing.  Mth  increased  Reynolds  number,  the  shear  layer  roll  up  is  triggered  at  a  shorter  streamwise 
location  resulting  in  a  near  field  which  spreads  more  rapidly  than  the  layer  in  the  low  Reynolds 
number  simulation. 

It  was  shown  that  the  streamwise  location  of  the  initial  vortex  pairing  scales  linearly  with  the 
log  of  the  amplitude  of  the  inlet  boundary  conditions.  Vortex  pairing  was  shown  to  be  the  domi¬ 
nant  mechanism  for  shear  layer  growth  even  up  to  1OOO0„.  The  effect  of  forcing  functions  on 
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tune-averaged  flow  properties  such  as  mean  velocity,  shear  layer  spread  rate,  and  turbulence  sta¬ 
tistics  were  investigated.  TTie  present  two-dunensional  simulations  captured  many  experimentally 
observed  phenomena.  Computed  spread  rates  agreed  very  well  with  those  from  experimental 
measurements.  Tune-averaged  properties  also  displayed  similar  features  observed  experimentally. 
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Comparison  of  v  -  velodty  spectrum  of  3*D  experiment^ 
that  of  present  simulations  using;  broad  mode  forcing :  -- 
and  discrete  (fundamentai)  mode  forcing : - . 


Vmlidty  contours  from  fundamental  mode  simulation  at 
ta2,7a8,if«  =  10*. 


Vortidty  contours  from  fundamental/first  subharmonic  mode  simulation 
at  fs  2.708, 


Vortidty  contours  firom  simulations  at  t  =  3.413,  L  =  20,  forced  with; 

(a)  broad  mode  only,  (b)  broad  and  fundamental  modes,  (c)  broad,  fundamen 
tal,  and  first  subharmonic  modes,  (d)  shear  layer  widths  from  (a)  •  (c), 

- :  case  (a) ; . :  case  (b) ; - :  case  {c\Re  =  10^. 


I 


FIG.  10.  Hme-averaged  properties  from  3  mode  simulation,  distance  from  inflow, 

X  =; - :  0.0 ;  +  ;  2.5  ;  5.1 ; . . . :  7.7 ;  <>  J 10-3 ; 

A  :  12.8 ;  □  :  15.4 ;  X  :  18.0,*e  =  10^. 
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FIG.  11.  Hme-averaged  properties  from  3<mode  simulation, 

- :  Oster  and  Wygnanski^  (Forced,  region  of  no  vortex  pairing), 

—  ■  •  •  —  :  Oster  and  Wygnanski^  (Forced,  region  of  vortex  pairing),  (c)  only 
Re  =  10^.  For  legend  see  Fig.  9. 


FIG.  12.  Time-averaged  properties  from  broad  mode  simulation,  (a)  •  (c)  Re  =  10^, 

(d)  •  (f)  =  10^,  distance  from  inflow,  x  =; - :  0.0 ;  -f  :  9.2 ;  ^  :  13.5 ; 

. . . :  17.8 ;  o  :  22.1 ;  A  :  26.4 ;  □  :  30.7 ;  X  :  35.0. 


FIG.  13.  Time-averaged  properties  from  broad  mode  simulation, 

- :  Oster  and  Wygnanski^  (Unforced,  naturally-developing  layer). 

For  legend  see  Fig.  11. 
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FIG.  14. 
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FIG.  15.  Streamwise  development  of  and  v^,  L  =  40, 

- ;  Jfe  =  10^ ; . ;  J?«  =  10^  (u^/^) 


t.o 


-•0  *20  20  60  100  140  180  220  260  300 


60  -20  20  60  100  140  100  220  260  300 


40  00  120  160  200  240 


220  *100  -140  -100  -60  -20  20  60  100 


00  -60  -20  20  00  ira  140 


-JOO  -260  -220  -100  -140  -100  -90  ~20  20  60 


-300  -260  -220  -100  -140  -100  -00  -20  20 


FIG.  16.  Rn  (O)  and R22  (A)  autocoiTdatioiis,  (a)  -  (e) Re  - 10^,  (f)  •  (j) ite  s  10^, 
at  streamwise  location,  = ;  (a),  (I)  2.9,  (b),  (g)  5.7,  (c),  (h)  14.3, 

(d),  (i)  25  J,  (e),  (j)  36.0, - :2fMf  degree  curve  fit  of  at  x  =0. 
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FIG.  17.  R22  two-p<fint  spatial  correlations,  (a)  Re  =  10^,  (b)  Re  =  10^, 

at  streamwise  location,X0  = ;  -4- :  2.9 ;  ^  :  5.7 ; . . . :  14.3 ;  O  •  >  A  * 
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